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ABSTRACT 

^■f-^ ' The problem of computing the pulse profiles from thermally emitting spots on the surface of a 

,__( ' neutro n star in general relativity is reconsidered. We show that it is possible to extend [Beloborodovl 

\ (|2002f ) approach to include (multiple) spots of finite size in different positions on the star surface. 

^SJ ■ Results for the pulse profiles are expressed by comparatively simple analytical formulas which involve 

^ ^ ' only elementary functions. 

d . Subject headings: relativity — stars: neutron — X-rays: stars 
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1. INTRODUCTION 



X-ray emission from isolated neutron stars (NSs), first detected in radio pulsars (PSRs), is now increasingly observed 
^ in other classes of sources, most of which are radio-silent or have radio properties much at variance with those of 
1^ . PSRs. They include the therrnally e mitting NSs (XDINSs; e.g. TuroUa 2009i) . the cen tral compact objects in supernova 
remn ants (CCOs; e.g lDe LTlcall2008[) . the magnetar ca ndidates (SGRs and AXPs; e.g. lMereghetti|[2008l : IRea fc Espositol 
1201 It ) and the rotating radio transients (RRaTs; e.g. lBurke-SDolaoill2012l ). 

With the exception of some PSRs, the X-ray emission of which is dominated by a non-thermal component of 
magnetospheric origin, the spectra of all other X-ray emitting, isolated NSs exhibit one (or more) thermal component 
Q ' which, most probably, originates at the star surface. Since pulsations are observed, thermal X-ray photons come 
^ either from a localized, heated region, like in SGRs/AXPs and PSRs, or from the entire cooling surface with an 
C/3 ■ inhomogeneous temperature distribution, like in XDINSs. In this respect the analysis of the observed pulse profiles in 
d ' different energy bands is bound to reveal mu ch on the surface thermal map of the NS, on the physical size and position 
of the emitting regions (s) in particular fe.g. lZane fc Turollal[200llAlbano et al.ll2010f) . 

The problem of modelling the pulse profiles of a rotating, thermally emitting NS, includi ng the effects of gravitational 
ray bending, is an old one and has been thoroughly addressed in the literature (e.g. iPechenick. Fta clas fc CohenI 

XTl'. [TMl [Pap[T995l : iPage fc Sarm entolll996l: iPsaltis. Ozel k. DePeol [2OOOI : lBeloborodovll200l iZane fcTiI^olla,200l . In 
ff^ ' particular, in their classic paper lPechenick. Ftaclas fc CohenI ()1983[ ) analvzed the emission from two antipodal, uniform. 



00 ' circular caps. Although their approach contains no inherent complexity, the treatment of photon propagation in 
\^ a Schwarzschild spacetime leads to elliptic integrals and requires numerical evaluation. In general, resorting to a 
^y-^ numerical approach is unavoidable every time a continuous surface temperature distribution, anisotropic emission 
and/or an arbitrary shape of the emission regions have to be accounted for. However. Beloborodov. (2002), by means 
, of a clever approximation, has shown that simple, analytical expressions can be derived for the pulse profiles in full 
I ■ general relativity (Schwarzschi ld spacetime) for po int-like spots. 

• ] In this paper we make use of lBeloborodovl ()2002|) approximate treatment to extend his analysis to the case of finite, 
uniform, circular spots. Our results are valid for an arbitrary number of spots, regardless of their size, temperature and 
mutual position on the star surface (e.g. two different, non-antipodal caps). Some more complex emission geometries 
(like a cap surrounded by a corona) can also be easily accommodated. The expression for the total observed flux is 
, analytical and this makes our approach both simple and fast for the evaluation and comparison of pulse proflles with 
observations. 

2. OBSERVED FLUX 

Let us consider a surface element dS on a neutron star of radius R and mass M and let us assume that the 
Schwarzschild solution correctly describes the spacetime outside the star (in the following Rs = 2GM/c^ is the 
Schwarzschild radius). Let us further introduce a spherical coordinate system, (r, 0, 0), centered on the star in such a 
way that the polar axis coincides with the line-of-sight (LOS; unit vector 1). The distance to the observer is Z? ^ i?. 

Because photon trajectories are not straight lines, the ray from dS = R^ sin 9d9d<j> which reaches the observer leaves 
the surface, with respect to the local normal, at an angle a 9 (see Figure [1]). The relation between a and 9 is given, 
implicitly, by the two equations 



dr 



6^ 7-2 \ 7. 



-1/2 



(1) 
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Figure 1. A schematic view of ray propagation. The angles 9, a and 9p are also shown. 

where b is the ray impact parameter (|Beloborodovll2003 ) . 
The (monochromatic) flux dF^, from dS detected by the observer is then 

,^ / Rs\ ^ ( dcosa\ dS 

where v is the photon frequency and the specific intensity, both measured by the static observer at r = i?. The 
total flux is obtained by integrating the previous expression over the visible part of the emitting region, Sy . If the 
emission is Planckian at the local (uniform) temperature T, = B^{T) and this results in 

In Newtonian gravity it is a = and the flux is simply proportional to the area of the visible emitting region projected 
i n the plane of the sky. 

IBeloborodovl (,,2002) found that a simple, approximate expression can be used to link a and 9, without the need to 
solve (numerically) eqs. ^ and ([2|), 

l-cosa = (l-cos6') (^1-^^ . (5) 

Eq. ([S|) is remarkably accurate and produces a fractional error < 3% for i? > 3Rs ■ Substituting cos a and d cos a /d cos 6 
into eq. Qi, one obtains 



2 



Sv 



R \ R 



The flux is, then, expressed by the sum of two contributions, the first proportional to the surface area and the second 
to the projected area of the visible part of the emitting region. The latter, apart for the factor (1 — R/Rs), is the 
analogue of the Newtonian expression, while the former is a purely relativistic correction. The problem of computing 
the flux, once the geometry is fixed, is therefore reduced to that of determining Sy and evaluating the two integrals 



/p= / cos 9 sine ded(j), I, ^ sm9d9d(j). (7) 
2.1. Single circular spot 

In order to proceed further, we consider first the simplest case, in which the emitting region is a circular cap of 
semi-aperture 9c with its center at (i?, 9o, 0). For the sake of simplicity, and also because this is the most common 
occurrence, we consider only the caseQ 9c < n/2. Moreover, we restrict to < 9o < tt, since the case n < 9o < 2tt is 
reduced to the previous one upon the substitution 9o ^ iir — 9o, given the axial symmetry around the LOS. 

The ^-integral in both Ip and Is is immediate. By denoting with </ib(0) the cap boundary (0 < 0b < tt), it is 



^ Despite this limitation many configurations of interest can be nevertheless treated (see i|2.2l l. 
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Figure 2. Geometry for a single spot; the visible part is in orange and the terminator is marked by the blue, dashed line. (Left) The 
spot is fully into view and 9o — 8c > 0. (Right) The spot is partially into view and do — Oc < 0; the darker part bounded by the dashed 
orange line is the region < 6 < l^o — 9c\ (see text). 



cos 9 sin e 4>b{e) de , Is ^2 



sine4>b{e) de , 



(8) 



where 6^^^, e^max 
The function ch 



are the limiting values of the co-latitude, which are discussed below. 

f, can be readily found noticing that a generic point on the cap boundary has coordinates {R , 9, 0f,). 
In a spherical coordinate system with the polar axis coincident with the cap axis (unit vector c), its coordinates are 
(i?,0 = 9c, the latter system is rotated by an angle 9o with respect to the former around an axis perpendicular 
to the 1-c plane. By exploiting the transformation between the (cartesian) coordinates in the two systems, one gets 



cos 9 — cos 9o cos 9c — sin 9o sin 9c cos $ . 



(9) 



sin 9 cos 0b = cos 9o sin 9c cos $ + sin 9o cos 9c ■ 
Solving the second for cos $ and substituting into the first one, one finally obtains 

cos 9c — cos 9o cos 9 

cos 06 = ■ a ■ a ' 

sm 9o sm 9 



(10) 



(11) 



cos (pb is actually related to the function h introduced bv lPechenick. Ftaclas fc CohenI ( 19831 ). It is immediate to verify 
that it is —1 < cos(/)f, < 1, and hence (j)b is defined, in the range 16*0 — 9c\ < 9 < min [6*0 + 9c, 27r — {9o + 9c)], since it 
must be, by definition, < 9 < tt. It is important to notice that the 0-range in which it is possible to define cos^b 
is not sufficient to cover the entire cap when the LOS intersects the cap itself: this occurs either for 9o — 9c < or 
9o + 9c > TT. Since it is < 0c ^ ■'r/2 by assumption and disregarding visibility, there is just one intersection, a,t 9 — 
or TT. In these cases the cap is fully covered only by adding the range < 9 < \9o — 9c\ or 27r — {9o + 9c) < 9 < tt, 
respectively. This accounts for the missing surface which is itself a circular cap, perpendicular to the LOS, where (f) 
spans the entire range [0 , 27r] (see fig. [2]). Accordingly, the definition of 4'b{9) can be continuously extended as 



cos t/r — cos t/n cos t 



sm Oo sm ( 



\9o -~9c\<9 < min [9o + 9c, In - {9o + 9c)] 

0<9<\9o-9c\ if9o-9c<0 

2tt - {9o + 9c) < 9 < tt a 9o + 9c>tt 



(12) 



to include all cases. 

At each visible point of the star surface it has to be cos a > 0, and the terminator lies precisely at a = tt/2. From 
eq. ([5]) it follows that the terminator co-latitude is given by 



cos9f = ( 1 - ^ 



(13) 



it is cos9f < and 9p > tt/2, as expected, since relativistic effects bring more than half the sphere into view. 
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In the case the cap is entirely visible, i.e. it does not intersect the terminator, its co-latitude is in the rangqj 
^min = rnax (0, Oo ~ 6c) < < 6'max — mill [27r — {Oo + 9c)^ Oq + Od- The presence of the terminator (at 7r/2 < 9f < tt) 
is easily accounted for by replacing 9min (^max), as defined above, with Op every time it is Op < Omin {Op < ^max)- 
Summarizing, it is 



min [max (0, 9o — ^c), dp] 



and 



t^max = mm (6*0 +Oc,ffF)- 

Turning to the evaluation of /p^s, both integrals become trivial for (jji, — n and yield 



TT cos sin 6* (i6' = 7r(sin^ f?2 — sin^ 6*1) , Is =2 tt sin 6 d9 — 277(00362 ~ cos 9 1) 



d9, 



for any pair of angles , ^2- In the opposite case, the two indefinite integrals 

cos 9c — cos 9o cos 9 



h 



cos a sm u arccos 



sm tin sm t 



and 



I2 = 2 sin 9 arccos 



cos 9c — cos 9o cos 9 



sm o„ sm f 



d9 



(14) 
(15) 

(16) 
(17) 
(18) 



have to be calculated. It turns out that I12 can be evaluated analytically in terms of elementary functions (see the 
Appendix for more details) 



Ii — sin^ 9 arccos 



cos 9c — cos 9o cos 9 
sin 9o sin 9 



and 



Bin Uq siii u 

— COS 9c \J — [cos 9 — cos{9o + 6*0)] [cos 



- sin^ 9c cos 9o arcsin 



cos 9 — cos 9o cos 9c 
sin 9o sin 9c 



9 - cos{9o - 9c)\ 



(19) 



COSf 



: —2 COS arccos 
+sign {9o + 9c 
—sign {9o — 9c) arcsin 



cos c/n cos f 



2 cos 9c arcsin 

sm tio sm f 

(cos 9o cos 6'c + 1) cos + sin'' 9o — cos^ 9c — cos 9o cos 



cosf 



cos tin cos t 



TT I arcsin 



sm Wo sm ( 

,2/1 „„„2, 



(1 + cos0)| sin^o sin^d 
(cos 9o cos 9c — 1) cos 9 + sin^ 9o — cos^ 9c — cos 9o cos ( 



(1 — cos 0)1 sin 6*0 sin 0c 



(20) 



where the arbitrary constant was set to zero. It is then Ip ^ = Ii^2{9max) ~ Ii,2{9min)- We note that Ji 2 take a simple 
form for 9 = 9o i: 9c- In particular, if the cap is fully into view (see Figure [2l left), it is Ip = i: cos 9 o s\t? 9 c and 
Is = 27r(l — cos 6'c), as it follows also from geometrical considerations. The complete form of eqs. ([T9|) and ([20| is 
actually required only when evaluating the integrals at 9p. 
The flux (eq. [6]) is finally written as 



= ( 1-^ ) B,iT)A,s{9c,9o), 



where we introduced the "effective" area 

^cff(0c,0o)-i?' 



i?5 
R 



1 



R 



(21) 



(22) 



2.2. Multiple spots and other geometries 

Having computed the flux seen by a distant observer for a single circular spot, it is straightforward to generalize the 
result to an arbitrary number of spots. We stress that this is possible because, using Beloborodov's approximation, the 
flux is proportional to the "effective" area of the cap, A^s, introduced in the previous section. Although we impose no 
restrictions on the parameters, the assumption that the spots do not intersect is understood. For the sake of simplicity, 
here we consider just two circular, uniform caps with semi-aperture 9c,i, and temperature Ti {i — 1, 2). Let us further 



* Note, however, that it cannot be 0„ 



■ and ^?max = 2n — (80 + Oc) because this would imply 6c > n/2. 
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Figure 3. (Left) Normalized flux as a function of phase for a single spot of semi-aperture 9c = 40° and different values of x E'nd 5; {right) 
same, but for 6c = 3°. 

assume that the spots are aligned, in the sense that their centers lie on the same meridian (the more general case 
of misaligned caps will be discussed in the next section), and let 6d,2 be the relative angular displacement of second 
one with respect to the first, the center of which is at do (of course it is 64^1 = 0). The spot centers are then at 
Oo,i — Oo + 0d,i and the total flux can be calculated by simply adding the two contributions 

Fj°^ = (1 - [MTi)AM,i,eo)+B,iT2)A,si0c,2,0o + 0d,2)] ■ (23) 

The case of a two-temperature cap, i.e. a cap at Ti surrounded by a circular corona at T2, is treated much in the 
same way by subtracting from the larger cap the contribution of the inner spot and adding the latter at the proper 
temperature 

F™^ = (1 " ^) {^-(^2) [AM0c,2,0o) - A,si9c,i,0o)]+ B,iTi)A,gi0cA,Oo)} . (24) 

With the aid of the previous expressions more configurations can be modelled. In particular, for a NS with a thermal 
map made of two (antipodal) caps at Ti while the rest of the surface is at T2, one obtains the flux by using twice eq. 
(|24p with 0c, 2 = 7''/2, the second time replacing 60 with 60 + tt, and summing the two contributions. The similar case 
of a single cap at Ti is handled by summing the flux given by eq. (p4)) with 9c 2 = 7r/2 and that of eq. ((2T|) . with 
= 7r/2, eo^0o + n. 

2.3. Pulse profiles 

In order to compute pulse profiles we consider first the case of a single spot. Let f be a unit vector parallel to 
the rotation axis and fl the star angular velocity, il = 2tt/P where P is the spin period. Observed periods in X-ray 
emitting INSs are in the range « 0.1-10 s, so the assumption of Schwarzschild spacetime previously introduced is fully 
justified. We also introduce the angles x, ^ between the LOS, the cap axis and the rotation axis, respectively, i.e. 
cos X = f • 1 and cos ^ = f • c. 

Since the cap co-rotates with the star, the vector c rotates around f , keeping ^ constant. This implies that 9o changes 
in time. Introducing the rotational phase 7 = fit + 70 (70 is an arbitrary initial phase), from simple geometrical 
considerations it follows that 



cos 9o = cos X cos ^ — sin x sin ^ cos 7 . (25) 

Eq. (1211) then provides the phase-resolved spectrum once the previous expression is used for cosSq. The pulse profile 
in a given energy band is immediately obtained integrating over frequencies. Since B^{T) dv = (jT^/tt (cr is the 
Stefan-Boltzmann constant), the pulse profile in the [t^i, 1^2] range is 

F{u^,u2)={l-^^ C{v,,v2)'^A,s[0c,9o). (26) 

where C(j^i, V2) — {ir/aT'^) f^^ B^{T) dv. Similar expressions hold for other geometries. 

Some examples are illustrated in Figure [3l |4] and [5] where the pulse profiles are shown for seven values of the angles 
X, S, in the range [0°, 90°], step 15°. Because eq. (|25|) is invariant by exchanging x and 5, only the 28 pulse profiles 
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Phase Phase 

Figure 5. Same as in fig. \3\ but for two difl'erent, non-antipodal caps on the same meridian {left) and with the second spot shifted in 
longitude by 45° (right). 

which are actuaUy diverse are shown. In all cases it is M = 1.4 M© and R = 15 km [R/Rs — 3.6), corresponding to 
9f — 112°, and the pulse profiles refer to the bolometric flux (i.e. C — 1), normalized to (i^max + ^niin)/2. Figure[3] 
shows the case of a single spot for 0c = 40° (left) and Oc = 3° (right), small enough to te treated as point-like (see 
The pulse profiles for two equal, antipodal (6*^ 2 — 180°) spots are illustrated in Figure HI again for 9^ = 30° (left) and 
Oc = 3° (right). Figure E] (left) refers to two non-antipodal (9^,2 = 120°), different (6*^1 30° ,6*^,2 ^ 45°, kTi ^ 0.4 
keV, kT2 — 1 keV) caps, while the right panel illustrates the same case, but with the second spot shifted in longitude 
by 45°. The latter is simply obtained by adding a constant phase-shift to 7 in eq. ((25)) when 9o refers to the second 
spot. 

3. DISCUSSION AND CONCLUSIONS 

In this investigation we revisited the problem of computing the pulse profiles from thermally emitting spots on the 
surface of a neutron star in general relativity. Our goal has been to develop a sim ple approach which can be readily used 
for a quantitative comparison of models with observations. iBeloborodovl ()2002l ). by means of a suitable approximation, 
was able to derive analytical expressions for the pulse profiles in full GR for point-like, equal, antipodal spots. However, 
if more realistic thermal configurations are to be accounted for, going beyond the point-like approximation becomes 
necessary. We have shown that it is possible to extend Beloborodov's approach to include (multiple) spots of finite 
size in different positions on the star surface. Results for the pulse profiles are expressed by comparatively simple 
analytical formulas which involve only elementary functions. 

A qualitative comparison between point-like and finite-size spots is provided by Figures |3] (single spot) and S] (two 
equal, antipodal spots); since 9c = Q produces a vanishing flux, 9c = 3° was used instead to simulate a (nearly) 
point- like spot (see below). Indeed, the pulse profiles in Figure |4] (right) appear very similar to those discussed by 
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Figure 6. Contour plot for the pulsed fraction as a function of the two angles x ^nd ^ (full lines) for two antipodal, equal spots with 
semi-aperture 40° {left) and 3° (right). The dashed contours show the analytical result and the rectangle bounded by the heavy lines (left 
panel) separates the five regions labelled I, II, III, IV (Beloborodov 2002, see text). 

iBeloborodovl pOOl see his fig. 4j3 and the four "types" he introduced (class I, II, III, IV) are clearly recognizable. 
This is better seen in Figure [6] (right), where the pulsed fraction, defined as PF = {F^ax — Fmin)l{Fmax + Fmin), is 
shown as a function of x and ^, together with Beloborodov's analytical result (his eq. [8]). The two sets of contours 
are nearly indistinguishable and the maximum pulsed fraction, equal to (i? — 2Rs)/{R + '2Rs) for point-like spots, is 
the same. 

As expected, for larger caps the pulse shape changes, the "plateau" disappears and the pulsed fraction decreases 
(figure m left). Now the constant PF contours are quite different with respect to those of a point-like spot, as clearly 
shown in the left panel of figure |6l The maximal pulsed fraction is ~ 30% lower than {R — 2Rs)/{R + 2Rs) for 
9c = 40°. In general, we find that the point-like approximation is reliable up to 9c ^ 5°. 

Larger caps can be treated either using the approach described here or resorting to methods based on general 
relativistic ray-tracing. We believe that the former offers a number of advantages, since it involves no numerical 
integration, and allows for a great flexibility, so that diverse thermal configurations of the NS surface can be modeled. 
An obvious limitation is that only purely blackbody (or at any rate isotropic) emission can be treated. Despite 
this simple model is often successfully used in fitting the (thermal components) of X-ray spectra, emission from the 
cooling surface of isolated neutron stars is expected to be more complicated, e.g. because the star is covered by 
an atmosphere, or because the emissivity is stron gly suppressed at energies below the electron plasma frequency if 
the surface layers are in condensed form (see e.g. iTurollal l2009l and references therein). Realistic emission models 
predict, to a different extent, an angular dependence of the emit ted intensity. While anisotropy is modest for non- 
magnetized atmospheric models fZavlin. Pavlov fc Shibanovlll996l) , it becomes substantial in magnetized atmospheres 
([Pavlov, Shibanov, Ventura & Zavlin 1994|), or in condensed surfaces (Turolla, Zane & Drake 2004). In general, it 
would be impossible to compute analytically the analogues of /1.2 (see ij2.ip for a non-isotropic radiation field. We 
point out, however, that i n cas e the intensity depends on the angle a only, i.e. = /jy(6') by using eq. (I5|), all the 
cosiderations presented in §2.11 still hold, although now numerical integration is required to obtain the pulse profiles. 
This, being the integral just over a single variable, 9, adds only a modest complication and the present method has 
still advantages with respect to fully numerical ray-tracing. Clearly this is not the case if Ii, depends on both angles, 
a and the associated azimuth, since double integrals should be evaluated. 

RT acknowledges financial support from an INAF 2011 PRIN grant. We are grateful to Paolo Esposito and Silvia 
Zane for their helpful comments on the manuscript. 
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APPENDIX 
INTEGRALS EVALUATION 
An integration by parts brings the first integral into the form 

1 r 



Ii = sin^ d arccos 



cos 9c — cos do cos 9 



sin 9o sin 9 



sin 9o sin 9c 



cos 9o — cos 6'c/i 
^^^^^^^^^^^^^= an 

2 



1- 



where fi = cos 9. By introducing z ~ (fi 

. \ COS 9c — COS 9o COS 

Ii = sm 9 arccos 



; — cos 9o cos 9c 
sin 9o sin 9c 



— cos 9o cos 9c)/ (sin 9o sin 9c 



sin 9o sin t 



sin^ 9c cos 9o [ , '^^ „ 



which, after some trivial manipulations, yields eq. (I19p . 
li is handled in a similar way. After integrating by parts, one gets 

„ „ r cos 9c — cos 9o cos 

h = —2 cost' arccos 



the previous expression becomes 

- cos 9c sin 9c sin 9o f '^^ 
J 



sin 9o sin 9 



/i(cos 9o — cos 9ciJi) 



sin 9o sin 9c 



Upon writing 




— cos 9o cos 9c ' 
sin 9o sin 9c 



u(cos 9o — cos 9cLl) „ if cos 9o + cos 9c cos 9o — cos ( 
= cos 9c-^ - ' 1 



eq. (jA3|) can be cast as 



1-^2 



1 + M 



(Al) 



(A2) 



(A3) 



(A4) 



'2 cos arccos 



cos 9r — cos 9r, cos ( 



sm Uo sm f 



2 cos 9r 



dz 



(cos 9o + cos 9c) / 

J (u + l) 



{fi + 1) ■\/— /i^ + 2 cos 9c cos 9o^i — cos^ 9c — cos^ 9o + l 

/d 
: 
— 1) ■\/— /i^ + 2 cos 9c cos 9ofJ. — cos^ 9c — cos^ 0o + 1 

The last two integrals in eq. (lASP are of the general type 

dx 1 . (6 — 2ap)x — bp + 2c 



{x + p)\/ ax^ + bx + c \Jbp— ap^ — c 



{x + p)\/b'^ — Aac 



(A5) 



(A6) 



(ap2 hp ^ c < 0,b^ — Aac > 0; e.g. iPrudnikov. Brvchkov fc Marichev! 119921 ) . No te, however, that the previous 
expression is valid only li x + p > 0. If a; + p < 0, as in the last integral in eq. (|A5P where /i — 1 < 0, a minus sign 
must be placed in front of the result. Making use of eq. (jA6[) and after some algebra, eq. I^U^ is finally recovered. 



